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Convection in the solar interior is thought to comprise structures 
on a spectrum of scales. This conclusion emerges from phenomen- 
ological studies and numerical simulations, though neither covers 
the proper range of dynamical parameters of solar convection. 
Here, we analyze observations of the wavefield in the solar photo- 
sphere using techniques of time-distance helioseismology to image 
flows in the solar interior. We downsample and synthesize 900 
billion wavefield observations to produce 3 billion cross-correla- 
tions, which we average and fit, measuring 5 million wave travel 
times. Using these travel times, we deduce the underlying flow 
systems and study their statistics to bound convective velocity 
magnitudes in the solar interior, as a function of depth and sphe- 
rical-harmonic degree €. Within the wavenumber band € < 60, con- 
vective velocities are 20-100 times weaker than current theoretical 
estimates. This constraint suggests the prevalence of a different 
paradigm of turbulence from that predicted by existing models, 
prompting the question: what mechanism transports the heat flux 
of a solar luminosity outwards? Advection is dominated by Coriolis 
forces for wavenumbers € < 60, with Rossby numbers smaller than 
approximately 1 0 -2 at r/R Q = 0.96, suggesting that the Sun may be 
a much faster rotator than previously thought, and that large-scale 
convection may be quasi-geostrophic. The fact that isorotation con- 
tours in the Sun are not coaligned with the axis of rotation sug- 
gests the presence of a latitudinal entropy gradient. 

imaging | thermal wind balance | Reynolds stresses | inverse problem 

T he thin photosphere of the Sun, where thermal transport is 
dominated by free-streaming radiation, shows a spectrum 
in which granulation and supergranulation are most prominent. 
Observed properties of granules, such as spatial scales, radiative 
intensity, and photospheric spectral-line formation are success- 
fully reproduced by numerical simulations (1, 2). In contrast, con- 
vection in the interior is not directly observable and likely 
governed by aspects more difficult to model, such as the integrity 
of descending plumes to diffusion and various instabilities (3). 
Further, solar convection is governed by extreme parameters 
(4) (Prandtl number approximately 10 -6 -10 -4 , Rayleigh number 
approximately 10 19 -10 24 , and Reynolds number approximately 
10 12 -10 16 ), which make fully resolved three-dimensional direct 
numerical simulations impossible for the foreseeable future. It 
is likewise difficult to reproduce them in laboratory experiments. 

Turning to phenomenology, mixing-length theory (MLT) is 
predicated on the assumption that parcels of fluid of a specified 
spatial and velocity scale transport heat over one length scale 
(termed the mixing length) and are then mixed in the new envir- 
onment. While this picture is simplistic (5), it has been successful 
in predicting aspects of solar structure as well as the dominant 
scale and magnitude of observed surface velocities. MLT posits 
a spatial convective scale that increases with depth (while velo- 
cities reduce) and coherent large scales of convection, termed 
giant cells. Simulations of anelastic global convection (6, 7, 8, 
9), more sophisticated than MLT, support the classical picture 
of a turbulent cascade. The anelastic spherical harmonic 
(ASH) simulations (6) solve the nonlinear compressible Na- 
vier-Stokes equations in the anelastic limit; i.e., where acoustic 
waves, which oscillate at very different time scales, are filtered 
out. Considerable effort has been spent in attempting surface 


(10) and interior detection (11, 12) of giant cells, but evidence 
supporting their existence has remained inconclusive. 

Results 

Here, we image the solar interior using time-distance helioseis- 
mology (11-13). Raw data in this analysis are line-of-sight photo- 
spheric Doppler velocities measured by the Helioseismic and 
Magnetic Imager (14) onboard the Solar Dynamics Observatory. 
Two-point correlations from temporal segments of length T of 
the observed Doppler wavefield velocities are formed and spa- 
tially averaged according to a deep-focusing geometry (15) 
(Figs. 1 and 2). We base the choice of T on estimates of convec- 
tive coherence time scales (6, 16, 17). These correlations are 
then fitted to a reference Gabor wavelet function (18) to obtain 
travel-time shifts 8t( 6, (p, T), where (6, (p) are colatitude and 
longitude on the observed solar disk. By construction, these time 
shifts are sensitive to different components of 3D vector flows; 
i.e., longitudinal, latitudinal, or radial, at specific depths of the 
solar interior ( r/R Q = 0.92, 0.96) and consequently, we denote 
individual flow components (longitudinal or latitudinal) by sca- 
lars. Each point (6, (p) on the travel-time map is constructed by 
correlating 600 pairs of points on opposing quadrants. A sample 
travel-time map is shown in Fig. 3. 

Waves are stochastically excited in the Sun, because of which 
the above correlation and travel-time measurements include 
components of incoherent wave noise, whose variance (19) di- 
minishes as T~ l . The variance of time shifts induced by convec- 
tive structures that retain their coherence over time scale T does 
not diminish as T~\ allowing us to distinguish them from noise. 
We may therefore describe the total travel-time variance o 2 ( T) = 
Tj 9 ,(P (<5t 2 (0, T)) as the sum of variances of signal S 2 and noise 

N 2 /T, assuming that S and N are statistically independent. 
Angled brackets denote ensemble averaging over measurements 
of 8t(0, cp, T ) from many independent segments of temporal 
length T. Given a coherence time T coh , we fit <j 2 (T) = S 2 + 
N 2 /T over T < T coh to obtain the integral upper limit S. The 
fraction of the observed travel-time variance that cannot be mod- 
eled as uncorrelated noise is therefore S 2 /o 2 (T coh ). For aver- 
aging lengths ^coh (=24 and 96 h) considered here, we find 
this signal to be small; i.e., S 2 « N 2 /T coh , which leads us to 
conclude that large-scale convective flows are weak in magnitude. 
Further, because surface supergranulation contributes to S, our 
estimates form an upper bound on ordered convective motions. 

Spatial scales on spherical surfaces are well characterized in 
spherical harmonic space: 

fafm(T)- [\in OdO [ 2 * d(P8z(0,(P, T)Y} m (0,<p), [1] 
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Fig. 1. Line-of-sight Doppler velocities are measured every 45 s at 4, 096 x 
4, 096 pixels on the solar photosphere by the Helioseismic and Magnetic Im- 
ager (background image). We cross correlate wavefield records of temporal 
length T at points on opposing quadrants (blue with blue or red with red). 
These "blue" and "red" correlations are separately averaged, respectively 
sensitive to longitudinal and latitudinal flow at (0,(p-,r/R Q =0.96), where 
(0, 4>) is the central point marked by a cross (see Fig. 2 for further illustration). 
The longitudinal measurement is sensitive to flows in that direction while the 
latitudinal measurement to flows along latitude. We create a travel-time 
maps 8t{6, <p, T ) by making this measurement about various central points 
(0, 4>) on the surface. Each travel time is obtained upon correlating the wave- 
field between 600 pairs of points distributed in azimuth. 


where Y lfm are spherical harmonics, (d,m) are spherical harmo- 
nic degree, and order, respectively, and 5z^ m {T) are spherical 
harmonic coefficients. Her e, we sp ecifically define the term 
“scale” to denote 2jiR q / + 1), which implies that small 
scales correspond to large € and vice versa. Note that a spatial 
ensemble of small convective structures such as a granules or in- 
tergranular lanes (e.g., as observed on the solar photosphere) can 
lead to a broad power spectrum that has both small scales and 
large scales. The power spectrum of an ensemble of small struc- 
tures, such as granulation patterns seen at the photosphere, 
leads to a broad distribution in *f, which we term here as scales. 
Travel-time shifts induced by a convective flow component 
v fm (r), are given in the single-scattering limit by = 

f Q r 2 dr 3£ f (r)v fm (r), where is the sensitivity of the measure- 
ment to that flow component. The variance of flow-induced time 
shifts at every scale is bounded by the variance of the signal in 
observed travel times; i.e., {8Sj m ) < S 2 /(j 2 (T coh )(Szj m (T coh )). 
To complete the analysis, we derive sensitivity kernels 3% f (r) that 
allow us to deduce flow components in the interior, given the as- 
sociated travel-time shifts (i.e., the inverse problem). 

The time-distance deep-focusing measurement (15) is cali- 
brated by linearly simulating waves propagating through spatially 
small flow perturbations, implanted at 500 randomly distributed 
(known) locations, on a spherical shell at a given interior depth 
(Fig. 4). This delta-populated flow system contains a full spec- 
trum; i.e., its power extends from small to large spherical harmo- 
nic degrees. The simulated data are then filtered both spatially 
and temporally in order to isolate waves that propagate to the 
specific depth of interest (termed phase-speed filtering). Travel 
times of these waves are then measured for focus depths the same 
as the depths of the features, and subsequently corrected for sto- 
chastic excitation noise (22). Note that these corrections may only 
be applied to simulated data — this is because we have full knowl- 
edge of the realization of sources that we put in. Longitudinal and 
radial flow perturbations are analyzed through separate simula- 
tions, giving us access to the full vector sensitivity of this measure- 



Fig. 2. The cross-correlation measurement geometry (upper box; arrow- 
heads — horizontal: longitude, and vertical: latitude) used to image the layer 
r/R Q = 0.96 (dot-dashed line). Doppler velocities of temporal length T mea- 
sured at the solar surface are cross correlated between point pairs at opposite 
ends of annular discs (colored red and blue); e.g., points on the innermost 
blue sector on the left are correlated with diagonally opposite points on > 
the outermost blue sector on the right. Six-hundred correlations are pre- o 

pared and averaged for each travel-time measurement. Travel times of waves o 

that propagate along paths in the direction of the horizontal and vertical £ 
arrows are primarily sensitive to longitudinal and latitudinal flows, v ^ and < 
v e , respectively. The focus point of these waves is at r//? 0 = 0.96 (lower 
box) and the measured travel-time shift 8x{0, (p, T ) is linearly related to 
the flow component v(r//? 0 = 0.96, 6 , (p) with a contribution from the inco- 
herent wave noise. We are thus able to map the flow field at specific depths 
v(r, 6, (f>) through appropriate measurements of 8r(0, (p, T). For the inversions 
here, we create travel-time maps of size 1 28 x 128 (see Fig. 3). For reference, 
we note that the base of the convection zone is located at r/R Q = 0.71 and 
the near-surface shear layer extends from r/R Q = 0.9 upwards. 


ment to flows. Travel-time maps from the simulations appear as a 
low-resolution version of the input perturbation map because of 
diffraction associated with finite wavelengths of acoustic waves 


Travel-time map (seconds) 



Fig. 3. A travel-time map consisting 1 6,384 travel-time measurements, span- 
ning a 60° x 60° region (at a resolution of 0.46875 ° per pixel) around the 
solar disk center, obtained by analyzing one day's worth of data taken by 
the Helioseismic and Magnetic Imager instrument (14) onboard the Solar 
Dynamics Observatory satellite. 3.2 billion wavefield measurements were 
analyzed to generate 1 0 million correlations, which were averaged and fitted 
to generate this travel-time map. This geometry and these particular wave 
times are so chosen as to be sensitive to flow systems in the solar interior. The 
spectrum of these travel times shows no interesting or anomalous peaks that 
meet the detection criteria (described subsequently). 
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Delta flows , size ~ X/4 




Fig. 4. Because wavelengths of helioseismic waves may be comparable to or larger than convective features through which they propagate, the ray approx- 
imation is inaccurate and finite-wavelength effects must be accounted for when modeling wave propagation in the Sun (20). In order to derive the 3D finite- 
frequency sensitivity function (kernel) associated with a travel-time measurement (21), we simulate waves propagating through a randomly scattered set of 
500 east-west-flow 'delta' functions, each of which is assigned a random sign so as not to induce a net flow signal (22) (upper box). We place these flow deltas in 
a latitudinal band of extent 120° centered about the equator, because the quality of observational data degrades outside of this region. We perform six 
simulations, with these deltas placed at a different depth in each instance, so as to sample the kernel at these radii. The bottom four boxes show slices 
at various radii of the sensitivity function for the measurement which attempts to resolve flows at r//? 0 = 0.96. Measurement sensitivity is seen to peak 
at the focus depth, a desirable quality, but contains near-surface lobes as well. Note that the volume integral of flows in the solar interior with this kernel 
function gives rise to the associated travel-time shift, which explains the units. 




excited in the Sun and in the simulations. The connection be- 
tween the two maps is primarily a function of spherical-harmonic 
degree €. To quantify the connection, both images are trans- 
formed and a linear regression is performed between coefficients 
of the two transforms at each € separately (see SI Appendix for 
details). The slope of this linear regression is the calibration fac- 
tor for degree €. 

We apply similar analyses to 27 d of data (one solar rotation) 
taken by the Helioseismic and Magnetic Imager from June-July 
2010. These images are tracked at the Carrington rotation rate, 
interpolated onto a fine latitude-longitude grid, smoothed with a 
Gaussian, and resampled at the same resolution as the simula- 
tions (0.46875 deg /pixel). The data are transformed to spherical 
harmonic space and temporal Fourier domain, phase-speed fil- 
tered (as described earlier) and transformed back to the real do- 
main. Cross-correlations and travel times are computed with the 
same programs as used on the simulations. Strips of 13 ° of long- 
itude and the full latitude range are extracted from each of the 
27 d results and combined into a synoptic map covering a solar 
rotation. The coefficients from the spherical harmonic transform 
of this map are converted, at each degree zf, by the calibration 
slope mentioned above, and a resultant flow spectrum is derived, 
as shown in Fig. 5. These form observational upper bounds on the 
magnitude of turbulent flows in the convection zone at the scales 
to which the measurements are sensitive. 

It is seen that constraints in Fig. 5 become poorer with greater 
imaging depth. This trend may be attributed to diffraction, which 
limits seismic spatial resolution to approximately a wavelength. In 
turn, the acoustic wavelength, proportional to sound speed, in- 
creases with depth. Because density also grows rapidly with depth, 
the velocity required to transport the heat flux of a solar lumin- 
osity decreases, a prediction echoed by all theories of solar con- 
vection. Thus we may reasonably conclude that the r/R Q = 0.96 
curve is also the upper bound for convective velocities at deeper 


layers in the convective zone (although the constraint at r/R Q = 
0.92 curve is weaker due to a coarser diffraction limit). Less re- 
strictive constraints obtained at depths r/R Q = 0.79,0.86 (whose 
quality is made worse by the poor signal-to-noise ratio) are not 
displayed here. 

Discussion 

Convective Transport. The spectral distribution of power due to an 
ensemble of convective structures, of spatial sizes small or large 
or both, will be broad. For example, it has been argued (10) that 
photospheric convection comprises only granules and supergra- 
nules, and that the power spectrum of an ensemble of these struc- 
tures would extend from the lowest to highest € . In other words, if 
granulation-related flow velocities were to be altered, the entire 
power spectrum would be affected. Thus the large scales which 
we image here (i.e., power for low € ), contain contributions from 
small and large structures alike, and represent, albeit in a com- 
plicated and incomplete manner, gross features of the transport 
mechanism. 

Our constraints show that for wavenumbers € < 60, flow ve- 
locities associated with solar convection ( r/R Q = 0.96) are sub- 
stantially smaller than current predictions. Alternately one may 
interpret the constraints as a statement that the temporal coher- 
ence of convective structures is substantially shorter than pre- 
dicted by current theories. Analysis of numerical simulations 
(6) of solar convection shows that a dominant fraction (approxi- 
mately 80%) of the heat transport is effected by the small scales, 
However, our observations show that the simulated velocities are 
substantially overestimated in the wavenumber band € < 60, pla- 
cing in question (based on the preceding argument) the entire 
predicted spectrum of convective flows and the conclusions de- 
rived thereof. We further state that we lack definitive knowledge 
on the energy-carrying scales in the convection zone. We may thus 
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Fig. 5. Observational bounds on flow magnitudes and the associated Rossby numbers. Boxes (A, B): solid curves with \ -o error bars (standard deviations) show g 

observational constraints on lateral flows averaged over m at radial depths, r/R Q = 0.92, 0.96; dot-dash lines are spectra from ASH convection simulations (6). k 

Colors differentiate between the focus depth of the measurement and coherence times. At a depth of r/R Q = 0.96, simulations of convection (6) show a 
coherence time of V coh = 24 hours (A) while MLT (16) gives T coh = 96 hours ( B ), the latter obtained by dividing the mixing length by the predicted velocity. 

Both MLT and simulations (23, 24) indicate a convective depth coherence over 1 .8 pressure scale heights, an input to our inversion. At r/R Q = 0.96, MLT predicts 
a 60 ms -1 , € = 61 convective flow and for r//? 0 = 0.92, an € = 33, 45 ms -1 flow [upon applying continuity considerations (23)]. (O shows upper bounds on 
Rossby number, Ro = U/{2£lL), L = 2nr/ ^/€{€ + 1), r = 0.92, 0.96 /?©. Interior convection appears to be strongly geostrophically balanced (i.e., rotationally 
dominated) on these scales. By construction, these measurements are sensitive to lateral flows i.e., longitudinal and latitudinal at these specific depths 
(r//? 0 = 0.92, 0.96) and consequently, we denote these flow components (longitudinal or latitudinal) by scalars. 


ask: how would this paradigm of turbulence affect extant theories 
of dynamo action? 

For example, consider the scenario discussed by Spruit in 
ref. 25, who envisaged very weak upflows, which, seeded at the 
base of the convection zone, grow to ever larger scales due to 
the decreasing density as they buoyantly rise. These flows are 
in mass balance with cool intergranular plumes which, formed 
at the photosphere, are squeezed ever more so as they plunge 
into the interior. Such a mechanism presupposes that these des- 
cending plumes fall nearly ballistically through the convection 
zone, almost as if a cold sleet, amid warm upwardly diffusing plas- 
ma. In this schema, individual structures associated with the 
transport process would elude detection because the upflows 
would be too weak and the downflows of too small a structural 
size. When viewed in terms of spherical harmonics, the associated 
velocities at large scales (i.e., low ), which contain contributions 
from both upflows and descending plumes, would also be small. 
Whatever mechanism may prevail, the stability of descending 
plumes at high Rayleigh and Reynolds numbers and very low 
Prandtl number is likely to play a central role (3, 25). 

Differential Rotation and Meridional Circulation. Differential rota- 
tion, a large-scale feature (€ ~ 2), is one individual global flow 
system and easily detected in our travel-time maps. Differential 
rotation is the only feature we “detect” within this wavenumber 
band. In other words, upon subtracting this € = 2 feature from 
the travel-time maps, the variance of the remnant falls roughly as 
T~\ where T is the temporal averaging length, suggesting the 
nonexistence of other structures at these scales. Consequently, 
we may assert that we do not see evidence for a “classical” inverse 
cascade that results in the production of a smooth distribution of 
scales. 


Current models of solar dynamo action posit that differential 
rotation drives the process of converting poloidal to toroidal flux. 
This process would result in a continuous loss of energy from the 
differentially rotating convective envelope and Reynolds’ stresses 
have long been thought of as a means to replenish and sustain the 
angular velocity gradient. The low Rossby numbers in our obser- 
vations indicate that turbulence is geostrophically arranged over 
wavenumbers € < 60 at the depth r/R Q = 0.96, further implying 
very weak Reynolds stresses. Because flow velocities are likely to 
become weaker with depth in the convection zone, the Rossby 
numbers will decrease correspondingly. At wavenumbers of 
€ ~ 2, the thermal wind balance equation describing geostrophic 
turbulence likely holds extremely well within most of the convec- 
tion zone: 


d£l _ C dS 
dz r 2 sin 6 36 ’ 


[ 2 ] 


where £2 0 is the mean solar rotation rate, Q is the differential ro- 
tation, z is the axis of rotation, 6 is the latitude, C is a constant, S 
is the azimuthally and temporally averaged entropy gradient. Dif- 
ferential rotation around € ~ 2 is helioseismically well con- 
strained; i.e., the left side of Eq. 2 is accurately known (e.g., 
ref. 26). The isorotation contours are not coaligned with the axis 
of rotation, yielding a nonzero left side of Eq. 2. Taylor-Proudman 
balance is broken and we may reasonably infer that the Sun does 
indeed possess a latitudinal entropy gradient, of a suitable form 
so as to sustain solar differential rotation (see e.g., refs. 27, 28). 

The inferred weakness of Reynolds stresses poses a problem to 
theories of meridional circulation, which rely on the former to 
effect angular momentum transport in order to sustain the latter. 
Very weak turbulent stresses would imply a correspondingly weak 
meridional circulation (e.g., ref. 29). 
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Let us define the forward and inverse spherical-harmonic transforms thus 



p27T PTl 

= / # do sin (9F/ m (<9,0) fo(0,0), 

Jo Jo 

L i 

(i) 

h(0,4>) 

£=0 m~-t 

(2) 


where (0, 0) are co-latitude and longitude on the sphere, Y^ m are spherical harmonics, t the degree and m 
order. 


Suppose we are interested in imaging a component p (longitudinal say) of sub-surface flows. Then we 
measure differential travel times Sr p (0,(f)) (Duvall et al. 1993; Hanasoge et al. 2010) that are specifically 
sensitive to this component. For flows whose magnitude is much smaller than the local sound speed, travel 
times Jt p ( 0,</>) are assumed to be linearly related to the underlying flow field through the following left 
convolution relation (Driscoll & Healy 1994), 

Sr p (n) = j r 2 dr J dll K p (TZ^n- r) • v(ft; r), (3) 

where 7 are elements in the Euclidean rotation group 50(3), and we map functions from 50(3) to 5 2 , 
the spherical surface. We denote the velocity field by v, and the kernel by K p , r is radius, superscript p 
refers to component p of the flow (i.e., longitudinal, latitudinal or radial) that we are attempting to image. 
Ideally, the measurement we make will only be sensitive to the desired flow component p but since other flow 
components may also affect the measurement (i.e., leakage), it therefore is more generally sensitive to the 
full vector flow field. Because we are studying the problem in terms of spherical convolutions, we effectively 
assume that sensitivity kernels are laterally invariant. 

First, we note that if we were to choose (for example) 


, a x ^ 5(r - r 0 ) 5(0 - 6 0 ) 5(cf> - fo) 

v<e, A r) = „„ e. 




where Vo is the amplitude of the delta function, then equation (3) becomes 

5t p ( 0, cf>) = Vo K p (7?. -1 (#o, 4>o)', fo) ■ e p , 


( 4 ) 

( 5 ) 
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i.e., the travel-time map is a rotated version (by an amount (0 q> </>o)) of the kernel itself (Duvall et al. 2006; 
Hanasoge et al. 2007) at radial location r = r o and where we define the delta function as 

/ r 2 drS(r ) = 1. (6) 

Jq 

In order to mimic wave excitation in the simulations as it is thought to occur in the Sun, we excite waves 
using a frequency band-limited stochastic driver (Hanasoge et al. 2006), implying that there is realization 
noise associated with measurements of travel times from a finite- length temporal calculation. To improve 
the signal-to- noise ratio of the retrieved travel times, we introduce 500 randomly distributed delta flows, 
compute the travel-time map and rotate the local map around each map and average them together to 
obtain the sensitivity kernel. 

The travel-time map in equation (3) is defined by three Euler angles. Consequently, a bijective mapping 
from 50(3) to S 2 is possible only when we consider a suitable subset of rotations. An example is when K is 
azimuthally symmetric (such as a Gaussian), a situation in which one of the rotations is rendered redundant. 
While K is not axisymmetric, we find, to first order, that the azimuthally symmetric (m = 0) component of 
the kernel is a dominant contributor, a conclusion that is further supported by the nearly linear curve the 
regression analysis (e.g., Figure 1) produces. Assuming that the kernel may be written as such, we invoke 
the convolution theorem (Driscoll & Healy 1994), which tells us that the convolution (3) of such a function 
in spherical harmonic space may be written as, 

ST !m = [ r2dr K ?o( r ) ’ Vlmir). (7) 

Jq 

In order to compute the kernel K^ 0 = K^, we transform the travel times and delta-function-populated 
velocity field into spherical-harmonic space (also see appendix ) and compute a linear regression between 
the two attendant sets of spherical harmonic coefficients at each t (e.g., see Figures 1 and 2). In a previous 
work (Hanasoge et al. 2010), we computed the sensitivity of the measurement used here to flows taken from 
a snapshot of ASH simulations. The ASH layer is radially coherent over a much broader depth range than 
the delta- flows used in the current calibration (effective radial width of 4.1 Mm). Empirically we find an 
approximate linear relationship between the measured response and assumed thickness of the convective 
layer (see Figure 3). 

The vector function = [ K J 1 ] is the sensitivity of the measurement aimed at imaging the p flow 
component of the underlying vector velocity field v^ m = [vl m (r)\. For example, if we were attempting 
to measure longitudinal flows, there may be leakage into the measurement from other flow components 
(latitudinal or radial). 

Now we consider the variance of the travel times, averaged over many ensembles (denoted by angular 
brackets) 

( ST em Sr em) = F f r2dr f r ' 2dr ' K T( r ) ( v lm( r ) v tm( r '))- ( 8 ) 

i,j J © J © 

We assume that flow velocities at two different depths are independent of each other, thereby allowing us to 
write 

(v\ m {r) vj m (r')) = ('■ v\ m (r)vl m {r )) 5{r - r'). 


(9) 




Fig. 1. — Linear regression analysis; spherical harmonic coefficients at t = 20, m E [—20,20] of the longitu- 
dinal velocity component ^ 20 , m on the x axis and St^q m on the y axis. The strong linear correlation between 
the two sets of coefficients, a relation that holds for all the coefficients, supports equation (7). The flows are 
tightly radially localized around t/Rq = 0.96. 


calibration vs. 1 for east-west flows, r/R=.96 



degree 1 


Fig. 2. — At each degree we perform a regression analysis akin to the one described in the caption 
to Figure 1 and derive one calibration constant. Here we plot this calibration constant as a function of 
spherical-harmonic degree t. It effectively describes the sensitivity of the measurement to longitudinal flows 
at r/R Q = 0.96 as a function of spatial scale of the flow. The error bars are at a 1 -a level (due to realization 
noise - introduced in the simulations in order to mimic the solar wavefield) - sensitivity is highest at the 
lowest degrees and drops with increasing wavenumber. The finite wave sensitivity up to £ ~ 65 is indicative 
of a full spectrum of flows (i.e., delta-populated flows) that extends to high degrees in £. 
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The variance now may be written as 

( Sr L ST em) = [ r 2dr K 7( r ) K e P ( r ) ( v em( r ) v em( r )) 

JQ 

+ E 

hj 

where the integral is split into two contributions, one due to the correlation between diagonal components 
of the kernels (the first term) and the other due to off-diagonal leakage between other components onto 
the desired measurement. Let us consider the case with p = (j) (for argument’s sake). Then, expanding 
equation (10) into its component terms, we have 

<[' Srfj 2 ) = I r 2 dr { [Kf] 2 ([vf„ J 2 > + [Kf] 2 <[<J 2 > + [iff] 2 <KJ 2 ) 

+ 2 Kf Kf (v r em v e im ) + 2 K? Kf {v\ m vfj + 2K++ Kf (vf m v ° tm ) } . (11) 

The first three terms on the right side are positive definite, and we find that (using data shown in Figure 4) 

Kf Kf, Kf * Kf r , Kf 4, Kf < [Kf 4 ] 2 . (12) 

Assuming that the cross-component Reynold’s stresses are comparable in magnitude to the kinetic energy 
terms, we have 

([Srfj 2 ) « [ r 2 dr [Kf 4 ] 2 ([vfj 2 ) + f r 2 dr [Kf] 2 ([- v e tm ] 2 ) + f r 2 dr [Kf] 2 (KJ 2 ), 

<[^/j 2 ) > [ r 2 dr [KfMvfj 2 ). (13) 

This argument holds for both lateral components, and we may write 

(Sr£ m Sr em ) > / r 2 dr K e (r) K e (r) (ve m (r)v em (r)), (14) 

Jq 

where we have dropped explicitly stating p and velocity components are written as scalars. Since this term is 
positive definite, we may weaken the inequality further and apply integration limits around a specific depth 
range D (say around r/R Q = 0.96) 

(StL)> [ r 2 dr K 2 (r) (v 2 m (r)), (15) 

J D 

We assume that the kernel is a separable function, i.e., Kj(r) = C^/(r), which is reasonable, since at depth, 
the kernel is radially localized and its horizontal structure is weakly dependent on depth. The term Ci is the 
calibration function that is obtained through the feature analysis. Defining 

Cim. = [ r 2 dr f(r) (v 2 tm (r)), (16) 

JD 

where Q m is the upper bound on the convective- velocity spectrum. We arrive at the desired inequality 

Ctm < (17) 

We show typical daily travel-time maps in Figure 5 and a synoptic solar rotation’s worth of daily maps in 
Figure 6. 


r 2 dr Kf (r) Kf (r) (v} m (r) vj m (r)) , 


( 10 ) 
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A. Details of the Method 


Since we use numerical methods (Hanasoge et al. 2006) to simulate wave propagation, we must use finite- 
sized “delta” functions. Consider the problem of imaging the convection zone in the region t/Rq G [0.75, 1.0]. 
Now, to compute the kernels and calibration curves to invert a given measurement to obtain the upper- 
bound flow spectrum, we need to determine the travel-time response to a delta flow at a series of depths. 
An individual “delta” function is defined as: 


(< 5 ) 

V cb = aeX P - 


9 • 2 

r sm 


OAcj) 2 


2 AO 2 


2a 2 


A r 2 \ 


(Al) 


Also, A (j) = cf) — <fi c , AO = 0 — 0 C , A r = r — r c , where (r c , 0 C , (j) c ) is the “location” of the anomaly and oy, cp 
are the radial and lateral widths of the flow deltas. We choose g\ = 3a r in our calculations. 


• The first set of flows sits between 0.9932 and 0.9999, r c = 0.9966, a r = 1.62, 

• The second sits between 0.9760 and 0.9932, r c = 0.9846, a r = 1.62 Mm, 

• The third sits between 0.9455 and 0.9760, r c = 0.9607, cr r = 1.74 Mm, 

• The fourth sits between 0.8995 and 0.9455, r c = 0.9205, a r = 2.96 Mm, 

• The fifth sits between 0.8364 and 0.8995, r c = 0.8659, a r = 2.96 Mm, 

• The sixth sits between 0.75 and 0.8364, r c = 0.7932, a r = 2.96 Mm. 

The choice of the radial and lateral widths of these features determines the calibration curve. The radial 
thickness of these features was varied in our experiment and we find an approximately linear relation between 
thickness and measured sensitivity (at least for the cases we considered). 

500 longitudinal (in one set of simulations) and radial (in a different set of simulations) flow deltas were 
randomly distributed at a given radial location. Each delta flow was also randomly assigned a sign and thus 
with this large sample of flow deltas, the net or total flow is nearly zero. The amplitude was chosen to be 
5% of the sound speed at r c , i.e., a = 0.05 c(r c ). 

The six data cubes listed above were generated. Then, to obtain the full 3D sensitivity kernel for each 
measurement geometry, the feature analysis was repeated on all six data cubes, giving us the sensitivity 
kernel at the depths listed above. So not only did we retrieve a calibration curve, but also extracted the full 
sensitivity kernels for all relevant measurement geometries. 
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Sensitivity to 1 km/s convective flow, r/R = 0.95 



Spherical harmonic degree, L 


Fig. 3. — Normalized measurement sensitivities of waves to a flows from an ASH simulation (Hanasoge 
et al. 2010) with a focus depth of r / R q = 0.95 and from the current study (red error bars; focus depth 
of v/Rq = 0.96). The ASH convective layer is radially much thicker than the delta- flows (4.1 Mm) used 
in the current calibration. We employ a simple normalization factor and divide the ASH sensitivity curve 
(Hanasoge et al. 2010) by the factor 22/4.1 = 5.4. For the range considered here, we empirically find an 
approximate linear relationship between the response and assumed convective layer thickness. 
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Fig. 4. — East-west travel time differences measured at the surface from simulations. Delta-function flows 
were introduced at various depths (denoted by r/R in red) and in the longitudinal and radial components. 
The time-distance geometry (‘antenna’) used has a focus depth of r/R® =0.96 and is primarily sensitive to 
east- west flows. All maps in columns C and D cover a longitude range of 330 degrees and latitude range of 101 
degrees about the equator. Image 7C shows the input at r/R = 0.96 displayed in the range ±3.5 km s -1 . 
There are 500 randomly placed small features in this image, each with the same amplitude and randomly 
picked signs. Cl - C6 are maps of travel-time shifts induced by longitudinal delta- flows introduced at specific 
depths. D1-D6 represent the same but with radial delta- flows at those depths. Columns C and D, lines 1-6 
are displayed in the range ±15 s. Cl-6 describe the sensitivity of the measurement designed specifically to 
isolate longitudinal flows to such flows whereas Dl-6 denote the undesired sensitivity to radial flow. In line 
with expectation, the largest responses are for signals that we are trying to extract, namely longitudinal 
flows at r/R = 0.96. To improve the signal to noise ratio, we average the local travel-time response about 
each of the (known) locations of the 500 individual features with the (known) signs taken into account. This 
improves the signal-to-noise ratio by a factor of V500 = 22. In column B(E) is shown the spatial averages 
about the 500 features for column C(D) with the scale ±5.4 s. Column A(F) shows the images of column B(E) 
scaled to a value appropriate to each image. For A1-A6, the ranges are ±0.7, ±0.7, ±5.4, ±5.4, ±0.7, ±0.7 
and for F1-F6, the ranges are ±3.5, ±1.3, ±0.5, ±0.4, ±0.3, ±0.3. The double feature in A4 presumably has 
to do with the opposite quadrants in the ‘antenna’. The relatively strong response in D1 is caused by the 
rays being nearly vertical near the surface and the sensitivity to vertical flows there. 




Fig. 5. — A travel-time map obtained by analyzing one day’s worth of data taken by the Helioseismic and 
Magnetic Imager instrument (Schou et al. 2011) onboard the Solar Dynamics Observatory satellite. The 
spectrum of these travel times shows no significant peaks which meet the detection criteria. 


Fig. 6. — A synoptic travel-time chart formed by extracting and sticking together 13.2° each from a solar 
rotation’s worth of daily maps (i.e., Figure 5). Like for Figure 5, the spectrum of these travel times shows 
no significant peaks which meet the detection criteria. 


